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Abstract 

We study the ABC model in the cyclic competition {A + B —* 2B, 
(N ; B + C -> 2C, C + A -> 2A) and the neutral drift (A + B -> 25 or 

2A, 5 + C -» 2C or 25, C + A -» 2A or 2(7) versions, with mu- 
, tations and migrations introduced into the model. When stochastic 

phenomena are taken into account, there are three distinct regimes 
in the model, (£) In the "fixation" regime, the first extinction time 
^C) • scales with the system size N and has an exponential distribution, 

\ with an exponent that depends on the mutation/migration probabil- 

ity per particle fx. (ii) In the "diversity" regime, the order parameter 
remains nonzero for very long times, and becomes zero only rarely, al- 
most never for large system sizes, (iii) In the critical regime, the first 
m ^ ' passage time for crossing the boundary (one of the populations becom- 

^ . ing zero) has a power law distribution with exponent —1. The critical 

mutation/migration probability scales with system size as N~ 1 . The 
transition corresponds to a crossover from diffusive behaviour to Gaus- 
sian fluctuations about a stable solution. The analytical results are 
checked against computer simulations of the model. 

1 Introduction 

Cyclic phenomena play a very important role in different classes of processes 
in nature, particularly in epidemiological and ecological systems. In the 
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epidemiological context, examples are diseases which do not leave us with 
permanent immunity PjJ |2] . In ecology, cases when three variants of a species 
compete with one-another in a cyclic fashion have been observed [HI EJ 
E]. Another system of interest are cyclic food webs. Recently, Kerr et al. 
obtained laboratory results for a system of three competing species that 
play rock-paper-scissors, when placed in a lattice-like spatial structure, and 
showed that there is agreement of the simulation results to the laboratory 
ones |7j. Its alternative, as far as evolutionary processes are concerned, is 
the famous Kimura- Weiss model of neutral genetic drift [HI EI - 

In a previous article JU] we have considered a non-spatial version of the 
ABC model in both the cyclic competition and neutral genetic drift versions, 
and studied the evolutionary time behaviour of such a model. The number of 
the A,B,C species oscillates with an amplitude that drifts with time, until 
one of the species (and then the second one) goes extinct, i.e. biodiversity is 
lost. The number of survivors vs. time plots show an exponential decay. It 
is very interesting to note that there is no difference in the extinction time 
scale for the ensemble of competing species and the case of neutral drift. This 
allows us to think of the neutral model as an "adiabatic approximation" of 
the cyclic system. Considering that the cyclic competition system would be 
a minimal model (three alleles) of Darwinian evolution picture, it seems that, 
by just looking at the results, we can not tell which one must have been the 
mechanism for the evolution! 

There is growing concern about the effects of habitat fragmentation in 
the survival of a species [TT]. Small, isolated forest fragments lose species-we 
have many examples in the case of fast deforestation worldwide |12j . Some 
species were lost within a few years of their isolation from the once-continuous 
forest ^3] ■ Oscillations in the number of population of one "species" are ob- 
served in the case of epidemics, as well as ecology. Throughout the past 
century we have observed changes in patterns of epidemics [TJj. For some 
diseases, the major transitions have been between regular cycles and irregular 
epidemics, and from regionally synchronized oscillations to complex, spatially 
incoherent epidemics. Sinervo et al. HH] observed spatio-temporal oscilla- 
tions in the number of lizards (male and female) that employ different mating 
strategies. 

Many population ecologists then conclude that a very important issue is 
the synchrony of population dynamics in different habitat patches [Td^ ITT| 
HHJ 120] ■ If the population of a certain species goes extinct in one patch 
while it still survives in other patches, then the hope is that what is known as 
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"rescue effect" can prevent global extinction [T7j. Otherwise, the population 
ecologists say, if there is synchrony of population dynamics between patches, 
the species is doomed to go extinct altogether. In this framework, there have 
been many proposed strategies for preservation, one of the most debated 
of which is the construction of "conservation corridors" that would make it 
possible for individuals to move along habitat patches |2H 122] • 

Considering that each habitat patch can be thought of as one copy of 
our non-spatial ensemble, it seemed reasonable to study the behaviour of the 
three-species cyclic system in the presence of mutations or migration. In the 
following section we introduce our model, and in the next ones, we present 
the results of our calculations and simulations. 



2 The Model 

Consider a non-spatial (well-stirred) system in which three species A, B, C 
are competing in a way described by the reactions: A + B — > 2B at a rate 
AB/N, B + C -> 2C, C + A -> 2 A at corresponding rates. Add to it the 
reactions: A — > B at rate fj,A, A — > C at rate fiA, B — > A at rate fiB, 
B — ► C at rate /X.B, C — > A at rate /xC, and C — ► at rate /xC, where /x is 
the probability of mutations per individual particle in unit time. 
The rate equations for this system will be: 
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+ fiB + fiC - 2fiA 

+ ^A + fiC - 2fiB (1) 
+ fiA + fiB - 2//C 

with A + B + C = N = const. If the neutral drift system is considered 
instead, i.e. when A + B — > 2A or 2B with equal probability (and similarly 
for the two other reactions), the rate equations will only contain the terms 
that depend on //, and the first two terms will be absent. 

In another version, mutations are replaced by migration into and out of 
the "island". In other words, the following migrations are happening: A 
out at rate 3/j.A, B out at rate 3fiB, C out at rate 3/iC, A, B, or C in at 
constant rate /x • N, where iV is the system size (to make it comparable to 
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the rate of leaving the island. This would correspond to the average number 
of individuals in other patches, which are the source of our "immigrants".) 
The rate equations will be the same as those for the system with mutations. 
The above equations have a fixed point at A = B = C = N/3, and it is a 
stable solution. 

However, the rate equations are just a "mean field" approximation; they 
only describe the behaviour of the average values of the individual popula- 
tions. In the real world, the system is subject to stochastic noise due to birth 
and death processes (intrinsic noise), which we take to be Poisson-distributed. 
The random nature of these processes need be taken into consideration, if 
we want to get the correct picture of the evolution of the system. For that 
we ought to write the master equation, and then try to somehow expand it, 
obtaining a Fokker-Planck equation, and eventually solve it. The two most 
commonly used expansions of the master equation are the Kramers-Moyal 
expansion, which is essentially a Taylor expansion in powers of system size 
N, and the f2-expansion of van Kampen j2H| which is an expansion in powers 
of \/N. The first method produces a nonlinear Fokker-Planck equation. The 
second method is quite systematic, and yields satisfactory results when the 
system has a single stable point, as the rate equations suggest is the case with 
the present model. In the absence of this stable solution, the task becomes 
very difficult [TUj . 

Using the "shift" operators notation: 



e A f(A,B,C) = f(A+l,B,C) 
e A 1 f(A,B,C) = f(A-l,B,C) 

the master equation for the cyclic competition system with mutations reads: 

dP(A,B,C,t) = { l_ [{ece -i _ 1)AC + {eAe - B i _ l)AB + {tBtc i _ 1)BC] + 

+^[{e A e~ B l + saSc 1 - 2) A + (e^ 1 + e B e A l - 2)B+ 

+ {e c e~ A l + eceB 1 - 2)C]}P(A, B, C, t) (2) 
while that for the neutral drift system with mutations: 

dP(A,B,C,t) _ 1 , , 

dt ~~ *-2N" +€a€ c -2)AC + 
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+ (W + egej 1 - 2)AB + (e^ 1 + t c t B x - 2)BC] + 
+fi[{e A e B 1 + e A e c l - 2) A + {e B e c l + e B e A l - 2)B+ 
+ (ecef + ece^ - 2)C)}P(A, B, C, t) 



(3) 



and that for the cyclic system with migrations: 



dP(A,B,C, t) 
di 



l)AC+(t A e B l -l)AB+(e B ec 



i 



1)BC]+ 



+M(e A - I) A + (e B - 1)B + (e c - 1)C+ 
+ j(e A 1 + e B 1 + ec 1 -3)}}P(A,B,C,t) 



(4) 



3 The 



Fixation" Regime 



In a previous work [IQ] we studied the cyclic competition and neutral ge- 
netic drift systems in absence of mutations. In the mean-field approximation 
(rate equations) the cyclic system has an infinity of neutrally stable solu- 
tions: any trajectory that conserves the integral H = ABC/N 3 is a stable 
trajectory of the system. On the other hand, the neutral drift system has 
an infinity of neutrally stable points: any state of the system is stable in the 
mean-field picture. The story is different, when stochastic birth- and-death 
processes are taken into account: then the individual population distribution 
drifts, until one of the species, and then another one, go extinct. In other 
words, there is fixation of the population to one of the varieties, and the 
product H = ABC/N 3 becomes zero. The survival probability versus time 
(in units of N) plot exhibits exponential decay with slope —3, which was 
verified by numerically solving the corresponding Fokker-Planck equation. 
This equation was obtained from a Kramers-Moyal expansion, since the in- 
finite multiplicity of stable solutions (points or limit cycles) makes the van 
Kampen method impractical to apply. The first extinction time scales with 
the total population size N. This behaviour is the same for both the cyclic 
competition and the neutral drift models, the second one behaving as an 
"adiabatic approximation" to the first one. For both models, the probability 
distribution function becomes uniform within a short time, i.e. any point 
inside the triangle that represents the phase space of the system is equally 
probable. 
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When /i is small (the meaning of small will become clear in section 5) 
the system with mutations/migrations is in the "fixation" regime, in which 
the results of our previous study [TU] are applicable qualitatively. In that 
regime, the stochastic processes push the system towards the (absorbing) 
boundary, in which at least one of the species has met extinction. This was 
verified by computer simulations of the master equations (0 El EI)- These 
simulations started with equal individual populations of A, B, C, (i.e. the 
centre). We generated times for the next possible reaction event with ex- 
ponential distribution as — tsirp ; where the rate of the cyclic/neutral or 
mutations/migrations as in the master equations above (0 HH EJ is substi- 
tuted. (Here rn is a random variable with uniform distribution in [0, 1]. This 
way we get Poisson distribution for the event times, i.e. really independent 
events [24 .) The reaction which occurs first is then picked and the system is 
updated accordingly. The process is repeated for a large number of events. 

We observed a "fixation" regime, in which one of the populations, and 
then the next one, go to zero quite fast. In this regime, the number of 
individual populations oscillates with an amplitude that drifts with time. 
The main reaction is then the cyclic one, and there are only occasional mu- 
tations/migrations, which are not frequent enough to prevent fixation. In 
other words, the system is experiencing "forces" due to fluctuations, and 
"forces" due to mutations (migrations). The fluctuations are pushing the 
system toward the boundary, while the mutations (migrations) push toward 
the centre. For small mutation/ migration rate the "fluctuation-generated 
forces" are stronger. Fig. Q shows the variation with time of the population 
of C's for a realisation of the system in the "extinction" regime. Occasionally, 
an individual of the competing species (the one next in the "food chain") is 
introduced in the system by mutations (or migrations), and then it either 
causes an occasional spike in the time series (like the one we see in Fig. ^ at 
25 < t < 30, or it "eats up" the old species completely and replaces it in the 
system (which is what happens between t = 35 and t = 40 in our plot). 

We investigated the behaviour of the probability distribution for the first 
crossing of the boundary, i.e. when one of the populations becomes zero 
for the first time, as the mutation/ migration rate increases. For that we 
ran ten thousand copies of the system for each value of the mutation rate. 
The survival probability (as defined above) was plotted vs. time. In the 
"fixation regime" the decay is still exponential, but the exponent varies with 
the mutation/migration rate, getting close to zero as the system approaches 
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criticality. Fig. |21 shows a plot of the exponent as a function of mutation 
probability \x for system size N = 150 and N = 300. 

4 The "Diversity" Regime 

If only mutations (migrations) are present, the system remains near the cen- 
tre point (A = B = C = N/3), and the order parameter remains consider- 
ably above zero; in other words, all three species are present in the system. 
One can occasionally observe temporary extinctions, but this happens very 
seldom, and for very long times; when the system size is very large, it al- 
most never happens. Since the boundary is not absorbing, occasional muta- 
tions/migrations will return the system to the state with maximal symmetry 
(biodiversity) where all three varieties coexist. When both the cyclic/neutral 
drift mechanism and mutations (migrations) are present, and the /i is above 
critical, the mutations manage to keep the system maximally disordered, 
since they are stronger than the fluctuations (which, as we saw, try to drive 
the system toward the boundary, i.e. fixation, and keep it there). Compared 
to the situation with no/low mutations (migrations), where the boundary 
is (practically) absorbing, and the final state of the system is a "pure" one 
(with only one of the species present), the mutation/migration rate acts then 
as a "temperature". 

Work has been done on the two allele almost neutral drift model with 
mutations [25J. The almost neutral model with mutations, preserving the 
total number of individuals, has only one degree of freedom, and allows one 
to derive an "effective potential" from the Fokker-Planck equation, obtained 
by a Kramers-Moyal expansion of the master equation. For small mutation 
probabilities, such that 2fiN <C 1, there is extinction of one species and 
fixation. The effective potential is symmetric around the centre (where both 
species are in equal numbers) and the branches of the effective potential 
are down. This allows for the system to quickly slip into the state where 
only one of the species is present. Otherwise, both species coexist forever 
in the high mutation regime, i.e. when 2fiN ^> 1. In that regime, the 
effective potential is symmetric around the centre point, but with branches 
upwards, which means that the centre point is a minimum potential point. 
The system then remains in the vicinity of that point for very long times. The 
effective potential "flips" from "branches up" to "branches down" at the point 
where 2fiN = 1. The transition is second-order, and critical behaviour is 
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observed. The system behaves similarly when migrations are present, instead 
of mutations. One aspect of migrations in a four-species system has been 
treated recently by Togashi and Kaneko |2o] . 

Since the rate equations have a stable solution, we employ the van Kam- 
pen expansion [23] . The idea of this expansion is to split the variables of the 
problem into a non-fluctuating part, and a fluctuating one, i.e. deal sepa- 
rately with the mean-field solutions and the fluctuations. In this approach, 
the numbers of the individual populations would be written as: 

A = Nfa + \fNx x 

B = N<J) 2 + VNx 2 (5) 
C = N(j) 3 + VNx 3 

Here the <f>i are the concentrations of A, B, and C species respectively (which 
only depend on time), and the the fluctuations (proportional to the 

square root of system size). The system size (total population) iV is con- 
served for the system with mutations, but not for that with migrations. This 
conservation rule will cause trouble in the case of the system with muta- 
tions, and will require special attention. Using the van Kampen Ansatz, the 
probability distribution P(A, B, C, t) is transformed into n({a;j}, t), and the 
following relations are true: 



U = N 3 / 2 P(N{4>i + VNx&t) 
dP _ 1 an 1 v d(pi Oil 



_l d_ 



1 d_ 1 d 2 

€l ~ ~JNdx~ + 2Ndx~^ + ''' 



Next we substitute everything into the master equation, leave only the 
term dU/dt in the left hand side, and group the right hand side terms ac- 
cording to powers of VN. The first term is of order iV 1 / 2 , and it must be 
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equal to zero, for an expansion in terms of N 1 ^ 2 to make sense. That term 
for e.g. the cyclic model is: 

^-^ dll d(pi 

]fa.\~tfjr + ^i^+i ~ 4>i<i>i+2 + /i(20i - (pi+i - 4>i+2)\ = 0, (7) 

which reproduces the rate equations in terms of the concentrations 0,, with 
steady state solution 0; = 1/3. Similarly, we get the rate equations for the 
other models. 

The terms of order N° give a linear Fokker-Planck equation of the form: 
<9n „ d , ^ 1 „ d 2 Tl . 

W = Eh^^(^n) + - 2 B*g^] (8) 

4.1 The System with Migrations 

We are going to solve the system with migrations first, since the absence of 
the conservation of total population TV makes this system easier to deal with. 
For simplicity, let us limit our attention to fluctuations around the steady 
state (pi = 1/3. The A-matrix for the system with migrations is: 

-3/i -1/3 1/3 

1/3 -3/x -1/3 

-1/3 1/3 -3/i 

and the B-matrix: 

2/t + 2/9 -1/9 -1/9 
-1/9 2/X + 2/9 -1/9 

-1/9 -1/9 2/X + 2/9 

These equations are linear, and the coefficients depend on time through 
(pi. This approximation is otherwise known as "linear noise approximation". 
The solution is known to be a Gaussian; the problem represents itself as an 
Ornstein-Uhlenbeck process. For our purposes, it suffices to determine the 
first and second moments of the fluctuations. Following van Kampen [2*3*] . 
we can multiply the Fokker-Planck equation by X{ and integrate by parts to 
get: 

(1/ *\ CCi ^ > X > a / „ \ 

dt =Y,A ij <x j > (9) 
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The eigenvalues of the A matrix are of the form — 3/i, — 3/x ± z/v3. The 
negativity of the eigenvalues guarantees the stability of the zero solutions to 
the first moments equations. Hence, the average of the fluctuations decays 
to zero and remains zero. The equations for the second moments can be 
obtained similarly: 

~ < X f^ 3 > =J2 A ik< X k%j > +Y1 4?& < XiX k > +Bij (10) 
k k 

By symmetry, in the steady state all the diagonal terms < x\ > are equal, 
as well as off-diagonal terms (correlations) < XiXj >. They depend on the 
migration probability \i alone. The steady state solutions are: 



2 9/x + l 
< xf >= — (11 



27/i 
1 



54// 



The diagonal terms coincide with the variance, since the < Xi >= 0. These 
fluctuations were measured "experimentally" , i.e. calculated from the results 
of the simulations. We simulated 1000 copies of the system (size N = 300) at 
different mutation rates, and calculated the mean and variance of the fluctu- 
ations, as well as correlations. The results of those simulations are shown in 
Fig-Elfor the system with migrations. They agree with the calculated values. 



4.2 The System with Mutations 

In the case of the system with mutations, the total population iV is a con- 
served quantity; this imposes restrictions on allowed fluctuations. We avoid 
this problem by excluding one of the variables. For that, we transform the 
problem into one with two variables, which can be done by putting: 
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This corresponds to transforming to Cartesian coordinates with the origin 
placed at the geometrical centre of the equilateral triangle, which constitutes 
the phase space of our system [TU] . 

The Fokker-Planck equation for the cyclic competition system becomes: 

The A-matrix has the form: 

-3/i -1/V3 \ 
1A/3 -3/i J 

with eigenvalues —3/i ± The mean value of the fluctuations decays 

to zero with an oscillating behaviour, and remains zero. The B-matrix is 
diagonal, with both diagonal elements equal to |(1 + 6/i). The equations for 
the second moments become: 



d < x 2 > 9 2 2 . 

—^ = -6/i<* >+-(! + 6/i) 

d < xy > 2 9 2 9 . . 

: = — = < x > -6u < xy > = < y > (14) 

dt V3 y V3 

d<y 2 > 2 „ 2 , 

with solutions that tend to 

2 2 1 + 6/i 

< x >=< y >= , < >= (15) 

as t — > oo. These solutions need to be transformed back in terms of < x\ > 
and < XiXj >, and the results are: 



1 + 6/i 



54/t 

For the neutral drift case, the Fokker-Planck equations becomes: 
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with the same steady state distributions as the cyclic competition model. 



5 The Transition Region 

The qualitative change in the behaviour of the system, when the value of 
the parameter \i is varied, speaks of the presence of a phase transition. The 
survival probability decays exponentially with time in the "fixation" regime, 
and the exponent goes down as the mutation probability increases, as shown 
in Fig. |21 This speaks of a symmetry-breaking transition as the mutation 
probability goes through the critical value. 

When the mutation/migration probabilities per particle approach zero, 
the leading term in the variances of the concentrations of individual popu- 
lations is of order 1/27N/J,, and it becomes of the same order of magnitude 
as the macroscopic concentrations, when 3fiN oc 1. This gives us the critical 
mutation/migration probability dependence on the system size N. The crit- 
ical /i was verified to be the same for both models with mutations present, 
and the critical \x for the system with migrations from the simulations is very 
close to that for the systems with mutations, which supports the calculations 
above. Fig. |U shows plots of critical \i vs. N, as well as a line of slope 1/N 
for comparison. To further verify this result, we performed a Kramers-Moyal 
expansion the models with mutations and migrations, in which \i was posed 
as Ho/N. Numerical solution of the resulting nonlinear Fokker-Planck equa- 
tions gave the equation for the eigenvalue (exponent of decay): A = _ 3, 
where A is the eigenvalue. At criticality, A = 0, which yields /Iq — 1/3 or 
/i = 1/3N [2*7] . This is in excellent agreement with the results obtained from 
the analysis above, as well as the computer simulations. It is worth noting 
that the neutral system essentially behaves like the cyclic one. 

The probability that the system has never crossed the boundary, (i.e. 
none of the populations has ever become zero,) when the system is in the 
transition regime, still decays to zero, but not exponentially any more. The 
transition is critical; however, transitions in nonequilibrium (steady-state) 
systems are different from the thermodynamic phase transitions. If the sys- 
tem size were infinite, one would introduce one infected individual (mutant) 
and measure the probability of survival of infection/mutation as time goes 
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to infinity. This would constitute the order parameter of our system in the 
thermodynamic limit. However, our systems are finite, and for finite N the 
quantity which exhibits critical behaviour is the first passage time for crossing 
the boundary 1 . 

In Fig. |3] we show plots of the first extinction times (with their cummu- 
lative probability in the y-axis) just below, just above, and at the critical 
fi (system size N = 210). (It is worth mentioning that we ran a statistics 
of 10000 copies of the system, but we are showing only 1 in 20 points in 
that plot, to prevent figures from becoming cluttered.) The power-law decay 
behaviour was used as a criterion for determining the critical point. The plot 
in Fig. shows clearly that the survival probability decays as a power law at 
the critical point. However, once the system is above the critical point, there 
is considerable probability for the system to have never "gone dead", even 
at very long times. The extinction times at the critical point scale with the 
system size. The power-law exponent is —1.03 ± 0.04 for cyclic competition 
model, —0.99 ± 0.03 for the neutral model with mutations; and —1.05 ± 0.06 
for the cyclic system with migrations. This exponent close to 1, found in sim- 
ulations is compatible with the one expected from branching processes |2*U] . 
as well as the one obtained for the two-species Kimura- Weiss model |2"5] . 

6 Conclusions 

We have considered an ABC model with cyclic competition/neutral drift 
and mutations (migrations) at a constant probability. The system exhibits 
a critical transition from a "fixation" regime to one in which biodiversity is 
preserved over long time. In the "fixation" regime, the number of the A, B, C 
species oscillates with an amplitude that drifts with time, until one of the 
species (and then the second one) goes extinct, i.e. the order parameter, de- 
fined as the product of ABC/N 3 goes to zero, and remains zero, except for 
occasional "bursts". In the "diversity" regime, the number of the A, B, C va- 
rieties fluctuates around the centre point, and there are rare extinctions, but 

^^Our model differs from the chemical reactions models: in them the onset of the cyclic 
behaviour is a Hopf bifurcation, in which a stable focus changes into a limit cycle |28| . In 
our model the limit cycle is absent. The chemical reactions models are dissipative even in 
the mean- field approximation, while our system has a centre in the mean field treatment. 
In those models the fluctuations become important only when the system is in the vicinity 
of the Hopf bifurcation, ours is entirely fluctuations-driven below the critical transition 
point. 



13 



the order parameter remains nonzero almost always. The survival probability 
decays exponentially below the transition point, but the exponent decreases 
as the mutation (migration) probability per particle increases, until it be- 
comes zero at the critical point. The critical mutation/migration probability 
depends on system size as N~ l , and the models have the same power-law 
exponent: -1. There is no difference in the behaviours of the neutral system 
and the cyclic system. Also, there is no qualitative difference between the 
system with mutations and that with migrations. 

These results address the concern about the effects of habitat fragmenta- 
tion in the survival of the species [11J. If the population of a certain species 
goes extinct in one patch (e.g. a herd, school, swarm) while it still survives 
in other patches, then it is hoped that the "rescue effect" can prevent global 
extinction [T71 [THJ CHI • Otherwise, the species is doomed to go extinct alto- 
gether. Our results suggest that the number of mutants/migrants necessary 
to preserve diversity is independent of system size, while in most realistic 
situations the number of mutations that happen would be proportional to 
system size. On the other hand, if one has in mind epidemiological systems, 
the important factor is to reduce migration probabilities below critical, which 
is exactly the purpose quarantines serve. Let's recall here the revival of SARS 
epidemics in Toronto, once the guard was let down, i.e. infected individuals 
were allowed to migrate from one community on to others. In that case, a 
large system size may be a disadvantage. It takes only one bad apple... 

In lattice models [3U1 small amounts of migration bring about "phase 
synchronization". Peak population abundances, however, are observed to 
be largely uncorrelated. Blasius et al. [T7j use a Langevin type system of 
equations to introduce noise in the three-species spatially structured model, 
and obtain complex chaotic travelling-wave structures. 

This work on three-species ecological or epidemiological systems relates 
to autocatalytic systems with a loop-like structure [321 12H]- Methods of anal- 
ysis employed in this paper can be expanded to the loop-like autocatalytic 
systems. We have work in progress regarding the four- or more-species au- 
tocatalytic loops 
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Figure 1: The time series for the number of C species in the "fixation" regime 
(here \i = 0.4 • 1CT 3 , and system size N = 300). 
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Figure 2: Dependence of the exponent of the decay of the survival probability 
on the mutation rate in the "extinction" regime for system size N — 150. 
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Figure 3: The variance and correlations of the fluctuations vs. migration 
probability per particle for system size N = 300, cyclic system. The same 
behaviour is observed for the system with mutations. 
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Figure 4: Variation of critical mutation and migration probability per particle 
with system size. 
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Figure 5: Cummulative survival probability vs. time (in units of N) just 
below, just above, and at the critical point, cyclic system with mutations, 
system size N = 210. 
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